AD-A218  454  FILE 


REPORT  DOCUMENTATION  PAGE 


Form  Appnvm 
OMB  NO.  0704-01 at 


rtoortmq  burden  lor  thu  coUwvo o  of  information  »t  nttmattd  to  avtr  aq*  I  "our  o*r  f  neons*,  mciuamq  tnt  tun*  tor  ravtrwinq  instruction*.  s«*r«imb  nmM  <jbU  sourcn 
gatfttnnq  $04  maintaining  tnt  dat*  090494,  and  comoifting  r«vt*wmg  tnt  <oH*cion  of  information.  Sand  commons*  rtgardinq  tmt  burdtn  rrtimat*  or  any  othtr  «*o*ct  of  this 
collection  of  information,  including  mggotttoot  for  reducing  tmt  ourotn,  to  Wawunqton  h*  adduarttn  $trvictt,  0 1  rector  at  t  for  information  Ooarauom  and  Attorn.  1215  JaHanon 
OavnMiqtmav.Swtt  1304,  Arlington,  v a  222Q2-*302.  and  to  tnt  Office  of  Management  and  Sudqtt.  Paoarwor*  A*ducttonAfO|*«t(0704^tU).w**Aington,  OC  2Q50J 


S^J  1.  AGENCY  USt  ONLY  (L«avt  Mink)  2.  REPORT  DATE  3,  REPORT  TYPE 

f>  ^  __  <  i_Final>Jlepa 

Q>  «.  title  and  suiTrru 

O  THE  EFFECTS  OF  CAUSTICS  IN  ACOUSTIC  INVERSE  SCATTERING 
EXPERIMENTS 

Z  AUTHORS) 

W.  W*  Symes 


3.  REPORT  TYPE  ANO  OATES  COVERED 

Final  Report,  1  Nov  88  thru  31  Oct  89 


5.  FUNDING  NUMBERS 


AFOSR-89-0056 
61102F  2304/A4 


|  7.  PERFORMING  ORGANIZATION  NAME(S)  ANO  AOORESS(ES) 

Rice  University 

Department  of  Mathematical  Sciences 
P.  0.  Box  1892 
Houston,  TX  77251-1892 


IffUN  mm  8.  PERFORMING  ORGANIZATION 

ftrvwit .  TR .  9«0°*T  3 


■ 


SPONSORING/ MONITORING  AGENCY  NAME(S)  ANO  AOORCSS(ES) 

AFOSR/NM 
Building  410 

Bolling  AFB,  DC  20332-6448 


SUPPLEMENTARY  NOTES 


10.  SPONSORING/ MONITORING 
AGENCY  REPORT  NUMRER 

AFOSR-89-0056 


M.  DISTRIBUTION /AVARAMUTY  STATEMENT 


12b.  DISTRIBUTION  COOf 


Approved  for  public  release  ; 
distribution  unlimited. 


13.  ABSTRACT  (Maximum  2QQ  words) 

i  Mo-tt  inversion  technique*  described  in  the  literature  rely  on  the 

validity  of  ray  tracing.  which  breaks  down  in  the  presence  of  caustics. 
. «.  . — s  The  linearized  acoustic  invert  problem  with  constant  reference  veloc* 

T  it y  is  analv zed  in  order  to  quantify  the  effects  of  a  caustic  in  a  probing 

S  f  1  B  S  *  wavefront  on  tin*  M*att**red  Ms»nal. 

.3.  A  When  the  sound  velomv  is  perturbed  bv  a  localized,  unidirec* 

ELECTE  BBA  tioual.  high  fr<vq«;*'iicv  mhomngen’itv,  the  Mirprising  result  obtained 
*  3*  IS  is  that  the  energv  m  the  mattered  field  is  spread  out  if  the  pertur- 

FEB  2  T  1990  f6  n  bat  son  is  located  *>n  the  c.vmk.  fins  >preading  of  energy  allows  the 

I  m  H  con>t ruction  of  an  i  Millatory  integral  representation  of  the  scattered 

DHmW  field,  which  has  «}*•*  Mine  form,  whether  or  not  an  incident  caustic  is 

present.  On  the  Hand,  a  sequence  of  localized  high  frequency 
sound  velocity  per»*:rhai*ons  constructed  such  that  the  size  of  the 
\J/V  scattered  signal  relative  r«>  the  *jze  of  the  mhomogeneity  becomes  arbi¬ 
trarily  large  as  » *  .*  *  ..f  the  perturbation  approaches  the  caustic. 

In  regions  •* « :e  jfp  « -unties,  a  general  iiiv’iu*  operator 

is  found  f«*r  .  vr.it  g  reference  velocities.  This  operator  Is 

>hown  to  He  ».(j  '»  operator  constructed  by  Boylkin 

I  ITmi 


•Terence  velocities.  This  operator  Is 
.'p*e  operator  constructed  by  Boylkin 


14.  SUBJECT  TERMS 

IS.  NUMBER  OF  FAGES 

78 

If.  PRICE  COOf 

17.  SECURITY  CLASSIFICATION 
OF  REPORT 

IS.  SECURITY  CLASSIFICATION 
OP  THIS  PAGE 

IB.  SECURITY  CLASSIFICATION 
OF  ABSTRACT 

20.  LIMITATION  OF  ABSTRACT 

UNCLASSIFIED 

.  UNCLASSIFIED 

UNCLASSIFIED 

SAR 

NSN  7540-01*280-5500 


Standard  Form  298  {Rav.  2*89) 

Prtvcnbtd  oy  AHU  Std. 


Final  Technical  Report 

AFOSR-89-0056 

The  Effect  of  Caustics  in  Acoustic  Inverse 
Scattering  Experiments 

W.W.  Symes 
February  2,  1990 

This  grant  supported  Dr.  Cheryl  Percell  in  her  final  year  of  graduate 
study  in  the  Department  of  Mathematical  Sciences  at  Rice  University.  Dr. 
Percell  completed  her  Ph.D.  work  under  my  direction  in  May  of  1989,  sub¬ 
mitting  a  dissertation  with  the  above  title  in  partial  fulfillment  of  the  degree 
requirements.  The  technical  accomplishments  of  her  research  are  summarized 
in  the  thesis  abstract,  reproduced  below.  A  copy  of  the  thesis  is  attached  to 
this  report. 


Abstract 

Most  inversion  techniques  described  in  the  literature  rely  on  the 
validity  of  ray  tracing,  which  breaks  down  in  the  presence  of  caustics. 
The  linearized  acoustic  inverse  problem  with  constant  reference  veloc¬ 
ity  is  analyzed  in  order  to  quantify  the  effects  of  a  caustic  in  a  probing 
wavefront  on  the  scattered  signal. 

When  the  sound  velocity  is  perturbed  by  a  localized,  unidirec¬ 
tional,  high  frequency  in  homogeneity,  the  surprising  result  obtained 
is  that  the  energy  in  the  scattered  field  is  spread  out  if  the  pertur¬ 
bation  is  located  on  the  caustic.  This  spreading  of  energy  allows  the 
construction  of  an  oscillatory  integral  representation  of  the  scattered 
field,  which  has  the  same  form,  whether  or  not  an  incident  caustic  is 
present.  On  the  other  hand,  a  sequence  of  localized  high  frequency 
sound  velocity  perturbations  is  constructed  such  that  the  size  of  the 
scattered  signal  relative  to  the  size  of  the  inhomogeneity  becomes  arbi¬ 
trarily  large  as  the  support  of  the  perturbation  approaches  the  caustic. 

In  regions  where  there  are  no  caustics,  a  general  inverse  operator 
is  found  for  smoothly  varying  reference  velocities.  This  operator  is 
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THE  EFFECT  OF  CAUSTICS  IN  ACOUSTIC 


INVERSE  SCATTERING  EXPERIMENTS 

CHERYL  BOSMAN  PERCELL 

ABSTRACT 

Most  inversion  techniques  described  in  the  literature  rely  on  the  validity  of  ray 
tracing,  which  breaks  down  in  the,  presence  of  caustics.  The  linearized  acoustic 
inverse  problem  with  constant  reference  velocity  is  analyzed  in  order  to  quantify, 
the  effects  of  a  caustic  in  a  probing  wavefront  on  the  scattered  signal. 

When  the  sound  velocity  is  perturbed  by  a  localized,  unidirectional,  high  fre¬ 
quency  inhomogeneity,  the  surprising  result  obtained  is  that  the  energy  in  the  scat¬ 
tered  field  is  spread  out  if  the  perturbation  is  located  on  the  caustic.  This  spreading 
of  energy  allows  the  construction  of  an  oscillatory  integral  representation  of  the  scat¬ 
tered  field,  which  has  the  same  form,  whether  or  not  an  incident  caustic  is  present. 
On  the  other  hand,  a  sequence  of  localized  high  frequency  sound  velocity  pertur¬ 
bations  is  constructed  such  that  the  size  of  the  scattered  signal  relative  to  the  size 
of  the  inhomogeneity  becomes  arbitrarily  large  as  the  support  of  the  perturbation 
approaches  the  caustic. 


In  regions  where  there  are  no  caustics,  a  general  inverse  operator  is  found  for 
smoothly  varying  reference  velocities.  This  operator  is  shown  to  be  equivalent  to 
an  inverse  operator  constructed  by  Beylkin  (1985). 
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CHAPTER  1 


Introduction 

Inverse  scattering  problems  are  central  in  many  data-processing  technologies, 
including  tomography,  radar  tracking,  seismology,  and  ultra-sonic  testing.  In  any 
inverse  scattering  problem,  if  the  medium  is  inhomogeneous,  then  wavefronts  gen¬ 
erally  will  become  singular.  The  variation  in  wave  velocity  causes  the  wavefronts  to 
fold  over  on  themselves,  as  in  Figure  1.1.  This  produces  an  imperfect  focusing  effect 
called  a  caustic.  In  this  case,  the  relationship  between  the  data  and  the  scatterer* 
is  not  well  understood  at  present.  In  order  to  begin  an  investigation  of  this  rela¬ 
tionship,  we  will  concentrate  our  efforts  on  the  linearized  acoustic  inverse  problem, 
which  is  directly  related  to  many  of  the  problems  listed  above. 

The  acoustic  wave  equation,  or  forward  problem  is: 

-st-t  uu  -  V2u  =  0  in  SRn  (1) 

c2(x) 

with  some  initial-boundary  data.  The  inverse  problem  can  be  stated  as  follows: 

Let  S(c]  =  u  on  (some  hypersurface  in  fRn)  x  (0,Tj.  Then 
given  Smeared’  find  c(x)  such  that  5(c]  =  ^measured  • 


1 


2 


to 


Figure  1.1. 


to 


(a)  An  inhomogeneous  medium  causes  a  wavefront  to  fold 
over  on  itself.  ( b )  The  caustic  is  the  envelope  of  the  singular 
wavefronts. 
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Linearization  of  the  inverse  problem  is  accomplished  by  a  perturbation  technique 
called  the  Born  approximation.  The  wavespeed  c(x)  is  assumed  to  be  a  small 
perturbation  about  a  smooth  background  velocity,  which  implies  a  perturbation  of 
the  acoustic  wave  field  u  about  a  reference  field 


c  =  Co  -f-  Sc 
u  =  Uo  +  6u  . 

Then  the  forward  perturbational  problem  can  be  expressed  as 

jrUo,t  ~  V2u0  =  0 

—  V2Su  =  fu0tt. 

The  linearized  inverse  problem  is  stated  as  follows: 


Given  Co,  let  6S  =  DS[cq\8c  =  6u  on  (some  hypersurface  in  3?n)  x 
[0,T].  Then  given  ^measured’  find  £c(x)  such  that  DS[cq]6c  = 
^measured1 


Most  inversion  techniques  in  practice  today  are  based  on  the  method  of  geomet¬ 
rical  optics.  The  behavior  of  a  wavefront  spreading  into  an  undisturbed  region  is 
obtained  by  propagating  discontinuities  along  characteristics,  or  rays.  The  theory  is 
valid  only  if  the  rays  fill  the  entire  region  of  interest  and  do  not  cross.  As  in  Figure 
1.2,  a  caustic  is  formed  by  the  envelope  of  crossing  rays,  making  the  geometrical 
optics  method,  and  hence  the  inversion  technique,  invalid. 


initial  wavefrout 


caustic 


Figure  1.2.  A  cusped  caustic  formed  by  a  convex  wavefront  in  a 
homogeneous  medium. 
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Even  though  a  caustic  is  defined  by  the  invalidity  of  the  geometrical  optics 
method,  it  is  a  real  physical  phenomenon  where  the  amplitude  becomes  very  large. 
In  fact,  Kulkarney  and  White  [2]  have  shown  that  in  a  medium  with  small  random 
perturbations  of  a  constant  wavespeed,  every  ray  passes  through  a  caustic  with 
probability  one  in  a  propagation  length  scale  of  order  0(a'~*),  where  the  standard 
deviation  of  propagation  speed  fluctuations  is  0(a). 

It  makes  sense  that  caustics  should  hold  some  information  about  the  medium 
and  perhaps  even  “illuminate”  it.  What  has  not  been  done  is  to  quantify  this 
information  and  take  it  into  account  when  solving  the  inverse  problem.  Perhaps 
Morawetz  [3,  p.  331]  said  it  best.  “There  is  as  much  information  about  the  sound 
speed  carried  in  a  caustic  as  in  a  smooth  ray  pattern.  The  problem  is  to  extract 
it.” 

Since  all  methods  available  today  for  solving  the  linearized  inverse  problem  ig¬ 
nore  the  presence  of  caustics,  it  is  natural  to  ask  what  effect  a  caustic  in  a  probing 
wavefront  has  on  the  scattered  signal,  ^u|surface*  Even  further,  we  want  to  know  if 
it  is  possible  to  represent  the  scattered  field  with  methods  that  are  available  today. 
Ultimately  the  goal  is  to  find  an  inverse  operator  for  that  representation. 

These  issues  are  addressed  in  the  following  manner:  In  Chapter  2,  we  review 
the  available  literature,  including  a  modification  to  the  geometrical  optics  method 
to  allow  solutions  near  caustics.  In  Chapter  3,  we  analyze  the  problem  of  a  convex 
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initial  wavefront  moving  into  a  region  with  constant  background  velocity.  The 
initial  wavefront  forms  a  caustic  and  the  perturbation  in  wavespeed  is  taken  to  be 
localized,  high  frequency  and  unidirectional  in  order  to  isolate  the  results.  The 
surprising  result  obtained  is  that  if  the  perturbation  is  located  on  the  caustic  in 
the  background  field,  the  scattered  signal  is  not  amplified;  it  is  in  fact  spread  out. 
The  problem  with  this  result  arises  in  taking  the  perturbation  to  be  more  general 
than  unidirectional.  In  Chapter  4  a  sequence  of  perturbations  is  constructed  so  that 
when  there  is  a  caustic  present  in  the  incident  field,  the  size  of  the  scattered  signal, 
6un ,  relative  to  that  of  the  scatterer,  8cn,  becomes  infinite  as  6c„  approaches  the 


caustic.  That  is 


MhsJE 

IIM2 


— ►  oo  as  n  — f  oo  . 


It  is  also 'shown  that  if  there  are  nr-  caustics  present,  then  jjf^Jjr  remains  bounded. 
In  Chapter  5  we  derive  an  inverse  operator  for  the  representation  of  the  scattered 
signal  obtained  in  Chapter  3,  and  in  Chapter  6  this  inverse  is  shown  to  be  related 
to  the  inverse  operator  constructed  by  Beylkin  in  [13], 


CHAPTER  2 


Background 


2.1  Geometrical  Optics 


Geometrical  optics  is  a  method  for  obtaining  the  behavior  of  a  wavefront  spread¬ 
ing  into  an  undisturbed  region  by  propagating  discontinuities  along  rays.  For  an 
excellent  discussion  of  this  method,  see  Whitham  [1]. 


The  wave  equation  is  solved  asymptotically  near  the  wavefront  by  means  of  an. 
ansatz  of  the  form 

u(x,  ~  £  aj(x,  O/jM*.  *)) 

;=i 

where  /j  =  /,■_  j. 

The  wavefront  is  given  by  <£(x,  t)  =  0,  and  it  is  assumed  that  u  is  identically 
zero  in  d(x,  t)  <  0.  For  us,  fj  has  either  the  form 


^m+.> 
{m+j) ! 


I  ®. 


<f>>  0 
4><  0 


for  a  wavefront  expansion,  where  Hq{4>)  is  the  Heaviside  function,  or  fi(<f>)  =  !»>» 
for  a  high  frequency  approximation. 
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Substitution  of  the  ansatz  into  the  wave  equation  (1)  and  setting  the  coefficients 
of  the  fj  equal  to  zero  gives  equations  for  the  phase  <f>  and  the  amplitudes  ay: 

|V^|2  =  p- <j> }  (eikonal  equation) 

\<f>tau  -  2 V<f>  •  Va!  +  -  V2<f>)a i  =  0  (1'*  transport  equation) 

c0  c 

The  eikonal  equation  is  solved  using  Hamilton-Jacobi  theory,  which  is  a  method 
of  characteristics  for  nonlinear  equations.  The  solution  is  obtained  by  solving  the 


following  system  of  O.D.E.s: 

p  =  V<f> 

T  =  <f>t 


fa  =P 

dt  _  r 
da  ~  ;*(X) 


S  =  -|Vco(x) 


In  this  case  <f>  is  constant  along  rays.  The  amplitude  ai(x,f)  is  also  given  by  the. 

solution  of  an  O.D.E.  along  the  rays: 

dax  dt  _  dx 

1 7  =  ati;*Va'Tc 
£  = 

These  equations  can  always  be  solved  as  long  as  the  rays  cover  the  entire  region 
of  interest  and  do  not  cross.  If  the  rays  tend  to  focus,  then  the  geometrical  optics 
approximation  incorrectly  predicts  that  as  the  ray  density  becomes  infinite,  so  does 
the  amplitude.  The  envelope  of  these  rays  is  called  a  caustic  (see  Figure  1.2). 

To  see  this,  multiply  the  transport  equation  by  at: 


T0  1 

a  ,  -r otai,  —  2Vd»  •  Vai  +  (-=■  <f>tt  —  V2^)aj 

[4  ^ 


=  0 


(3) 
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This  gives  the  space-time  divergence 


f  (5* '«■!)  +  V  •  (-a’V«  =  0 


Consider  a  slim  tube  formed  by  rays  emanating  from  the  initial  wavefront  in 
some  disk  Do  centered  at  a  point  Xo  and  going  to  the  wavefront  at  some  later  time 
t.  Let  Dt  be  the  slice  of  rays  at  time  t.  (See  Figure  2.1.)  Integrating  (3)  over  the 
volume  Rt  formed  by  the  rays  bounded  by  Do  and  Dt  and  using  the  divergence 


theorem,  we  get 


[  (4&,-a?V^).nds  =  0 

JdRt  Cq 


On  Do  the  normal  is  n  =  (—1,0);  on  £><,  n  =  (1,0);  and  on  the  sides  n  is  perpen¬ 
dicular  to  the  rays,  which  are  parallel  to  the  vector  ( \(f>t ,  —  V^).  Thus 


/  -  [  4<hds  =  0  . 

Jd,  c$  Jdo  c$ 


Suppose  x(X, t)  maps  a  point  on  the  initial  wavefront  to  a  point  x  at  time  t.  Let 
J(X,  t)  =  |V,x|  be  the  Jacobian  of  this  map.  Then 


/.[ 


«i(x(X, OvO  1  /  /V  r/Y  0)  1  /y  s}q _ n 

cg(x(X,0)  ^  cg(X)  -°- 


Since  this  holds  for  any  disk  Do  on  the  initial  wavefront,  then  as  long  as  the  integrand 
is  continuous, 


«.(X(X,0.|)  =  °.(x.0)a,-(^’i-))J-t(x,t) . 
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When  x(X,  t )  becomes  singular,  J  approaches  zero  and  (5)  is  no  longer  valid. 
But  expanding  (4)  to  first  order  about  Xo,  we  get 

c5(x„)AD'-4(Xo)AD°’ 

where  A D%  and  AZ?o  are  the  infinitesimal  areas  of  the  disks  Dt  and  Do.  This  is 
actually  a  statement  that  energy  flux  along  the  ray  tube  is  constant.  It  is  clear  that 
when  A Dt  approaches  zero,  a](x)  approaches  infinity.  The  ray  tube  collapses,  and 
the  theory  is  no  longer  valid. 


2.2  Caustics 

Although  geometrical  optics  fails  at  a  caustic,  it  can  be  altered  to  produce  valid 
results.  Expanding  on  the  example  of  Keller  [4],  Ludwig  [5]  has  constructed  uniform 
asymptotic  expansions  for  solutions  of  linear  hyperbolic  equations  at  a  caustic.  For 
the  simple  case  of  a  smooth  caustic  and  the  reduced  wave  equation 

V2u  +  (J02U  SS  0, 

the  ansatz  is  changed  to  a  linear  combination  of  the  Airy  function  and  its  derivative. 
The  structure  is  such  that  on  one  side  of  the  caustic  the  original  oscillatory  expansion 
is  given  by  geometrical  optics  (illuminated  side),  and  on  the  other  side  the  solution  is 
exponentially  decaying  (dark  side).  The  transition  between  the  two  sides  is  smooth 
with  asymptotically  larger  amplitude  along  the  caustic. 
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We  show  here  the  motivation  for  the  construction  since  it  will  also  motivate 
seme  of  our  own  analysis.  An  appropriate  ansatz  for  the  leading  order  asymptotic 
solution  of  the  reduced  wave  equation  is 

u(x)  =  a(x)e,w^  , 

where  u>  is  large. 

The  eikonal  and  transport  equations  are  respectively 

I W  =  1 

2V^  •  Va  +  a W  =  0  . 

4 

But  at  a  caustic,  the  local  ray  density  goes  to  zero,  the  amplitude  goes  to  infinity, 
and  the  geometrical  optics  solution  becomes  invalid. 

If,  however,  the  initial  wavefront  were  decomposed  into  plane  waves,  each  plane 
wave  could  propagate  into  the  region  without  forming  a  caustic  since  the  wavespeed 
is  identically  one.  This  suggests  an  ansatz  of  the  form 

u(x)  =  J  a{x,9)eik*Md9,  (6) 

where  <f>(x,  9)  =  constant  describes  a  plane. 

Ordinarily,  the  method  of  stationary  phase  can  be  used  to  get  ar„  asymptotic 
representation  for  -u.  If  the  stationary  points  defined  by 
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are  given  by  0j,  then 


u(x)  ~  £ 


2tt 


e±,'«a(x,^(x))e,'h'^) 


=  0  (u>~* ) 


as  long  as  ?W(x,0;)  £  0.  If  fafafy)  =  Mx>*j)  =  0  and  <f>999(x,$j)  -t  0,  then 
similar  analysis  yields  u(x)  ~  0(w"3  ).  If  two  stationary  points  coalesce,  then  <f>99 
approaches  zero,  and  in  the  region  between  <f>99  =  0  and  <f>99  ^  0,  the  method 
of  stationary  phase  breaks  down.  This  transition  zone  corresponds  to  the  region 
near  a  caustic.  Notice  in  Figure  2.2  that  each  stationary  point  corresponds  to  a 
contribution  from  a  ray.  Each  point  to  the  right  of  the  caustic  is  the  intersection 
of  two  rays.  Along  the  caustic,  each  point  corresponds  to  only  one  ray.  No  rays 
penetrate  to  the  left  of  the  caustic.  This  region  is  called  a  shadow  zone. 

The  simplest  function  modelling  this  phenomenon  is  given  by  Chester,  Friedman, 
and  Ursell  [6].  They  showed  that  £(x,  9),  j(x),  and  p(x)  can  be  found  so  that 

<£(x,  &)  =  7  +  p{-!j{3  , 

provided  that  <f>  is  analytic.  With  this  change  of  variables,  the  integral  (6)  becomes 

The  stationary  points  are  given  by 

rM-\e)^-e=o. 


u(x)-W«(x.«(fl) 
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As  -*  0,  the  second  derivative  also  goes  to  zero,  but  the  third  derivative  is 
constant. 


By  expanding  a(x, #(£));§  >n  a  Taylor  series 


a(xi0)^  =  ao  +  ai£  +  (£2  -  p)C?(£) 


and  integrating  by  parts  one  can  show  that 


u(x)  =  27re'Wx) 


ao(x) 


ir-Ai  p(x))  +  (-w>  p(x)) 

;3  x  '  i(jjs  x  ' 


+  0 


U1 


where 

A[-t)  =  T  ei[t^e)d£ 

2tt  J— oo 

is  an  Airy  function.  This  provides  a  new  ansatz 


ti(x)  =  eiun,{x) 


a0(x).4j  (-<*;»  p(x))  +  Q-f->4f/  (u>?  p(x )) 
v  7  v  7 


(7) 


Substituting  (7)  into  V2u  +  cj2u  =  0  and  collecting  like  powers  of  u>  gives  equations 
for  7(x),  p(x),  ao(x),  and  at(x)  corresponding  to  eikonal  and  transport  equations.  It 
can  be  shown  that  depending  on  the  sign  of  p(x),  the  eikonal  equation  is  hyperbolic 
(illuminated  region),  elliptic  (shadow  region),  and  parabolic  (caustic).  (See  Figure 
2.3.) 


In  the  illuminated  region,  the  original  geometrical  optics  approximation  is  still 
valid.  Near  the  caustic,  the  new  Airy  function  representation  is  appropriate.  If  the 


parabolic  illuminated 

caustic 


initial  wavefront 


Figure  2.3.  Depending  on  the  sign  of  p(x),  the  eikonal  equation  is 
hyperbolic,  elliptic,  or  parabolic. 
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caustic  is  analytic  and  a!  is  analytic  on  the  caustic,  then  7,  do,  dj  and  p  can  be 
continued  into  the  shadow  region,  and  the  solution  can  be  interpreted  in  terms  of 
complex  rays. 

There  are  similar  results  for  cusped  caustics,  linear  systems,  and  progressing 
wave  expansions  using  generalized  Airy  functions  or  solutions  to  the  Tricomi  equa¬ 
tion 

dx2f(x,y)  =  xdy2f(x,y) . 

Kravstov  [7,8]  has  also  developed  uniform  asymptotic  solutions  of  hyperbolic  equa¬ 
tions,  in  particular,  Maxwell’s  equations.  Stickler  et  al.  [9]  have  applied  Ludwig’s 
progressing  wave  ansatz  to  a  point  source  problem. 

2.3  Inverse  Problems 

There  are  many  methods  available  for  solving  linearized  inverse  problems.  Clayton 
and  Stolt  [10]  use  the  VVKBJ  method  to  asymptotically  obtain  the  Green’s  operators 
for  the  reference  field.  Then  a  comparison  of  the  Green’s  operators  to  the  kernel  of  a 
Fourier  transform  facilitates  the  relation  of  the  data  field  to  the  medium  parameters. 
Of  course,  the  WKBJ  method  is  invalid  in  the  vicinity  of  turning  points,  which 
correspond  to  points  on  a  caustic.  Cohen  and  Bleistein  [11]  use  a  family  of  specific 
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reference  fields,  that  they  call  probes,  to  solve  an  integral  equation  of  the  form 

J  u0va(r)dV  =  0 

for  the  medium  parameter  a(r).  Stacking  and  migration  methods  are  also  based 
on  ray  theory.  An  illuminating  exposition  on  the  connection  between  geometrical 
optics  and  migration  methods  is  found  in  a  paper  by  Lailly  [12]. 

Beylkin  [13]  gives  a  rigorous  derivation  of  migration  which  he  defines  as  a  dis¬ 
continuity  imaging  technique.  He  linearizes  the  inverse  problem  for  the  Helmholtz 
equation  with  a  point  source.  That  is, 

[V2  +  &2rc2(x)]v  =  S(x  —  7/) 

n2(x)  =  n2(x)  +  /(x) 

V  =  v'n  +  v3C 


gives 

[V2  +  k7nl]vin  =  6(x  —  rj) 

[V2  +  k7n20\vsc  =  -k7f(x)v'n  . 

Geometrical  optics  is  used  to  find  the  Green’s  function  for  v'n  and  for  vsc: 

i"  ~  a(k)Ain{x,rj)eiki'f'{x’,,) 
vout  ~  a(k)Aout(xi1^)eiK*,>'“(x,®  . 
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Then  the  scattered  field  is  represented  as 

»“(*,«> >l)  =  Wf(k,l,V)  =  (-it)”-1  f ' 

J  X 

where  (f>  =  4>in  +  <f>out ,  a  =  AmA0U*,  and  77  and  £  represent  source  and  receiver 
positions. 

The  solution  is  shown  to  be  related  to  the  Fourier  transform  of  a  generalized 
Radon  transform  (GRT) 

where 

=  J  f(x)a(x,{,q)6(t-t(x,(,tj))dx,  for  f>0 

77)  =  0,  for  t  <  0  . 

This  is  a  generalized  Radon  transform  since  {<f>  =  0}  can  be  any  hypersurface,  not 
just  a  plane. 

The  inverse  is  constructed  from  the  dual  to  the  GRT,  which  is  called  a  generalized 
backprojection  operator.  The  operator 

f/(y)  -  BTF*WHy)  =  ^  J~  lxJxem*MMMx,y,(,ri) 

xf(x)dKd£kn-ldk 

is  shown  to  be  asymptotic  to  the  identity  operator.  The  operator 
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is  applied  since  v,c  is  needed  in  the  time  domain.  Also, 


Since  the  Taylor  series  for  $  is 

$  =  V!/^(y,^^)-(x-y)  +  C?(|x-y|2) , 

the  following  change  of  variables  makes  F f(y)  look  almost  like  the  identity  operating 
on  f{y ): 


p  =  kvv<f>(y,S,n) 

dp  =  kn-lh(y,Z)d(dk  . 


Then 

where  nn( y )  is  the  part  of  the  region  where  no  caustics  are  present.  Of  course,  it  is 
critical  to  this  method  that 


h{y,0  =  det 


<?y,  <f>y2  .  .  .  <pyn 

Pyi<i  i  •••  OyuZi 


(8) 


^yi4n-l  1  '  •  *  ^.Vn5n-l 

be  non-zero.  This  is  equivalent  to  the  condition  that  there  be  no  caustics  present. 
In  fact,  Beylkin  uses  a  cutolf  function  to  ensure  that  h(y,£)  ^  0.  This  is  why  fl„(y) 


CHAPTER  3 


Plane  Wave  Velocity  Perturbation 

The  first  step  in  solving  an  inverse  problem  is  making  sure  we  can  solve  the 
forward  problem  associated  with  it.  A  potential  problem  is  that  in  the  right  hand 
side  of  the  perturbational  equation  is  the  product  of  Sc  and  u0tt.  One  might  expect 
that  a  perturbation  located  on  a  caustic  in  the  reference  field  would  be  illuminated. 
That  is,  it  would  have  an  effect  on  the  scattered  field.  It  is  not  clear  that  the  trace 
of  the  scattered  field  can  be  represented  using  geometric  optics  or  the  methods  of 
Ludwig. 

In  this  chapter  we  analyze  an  appropriate  problem  in  order  to  investigate  the 
effect  on  the  scattered  signal  of  having  a  caustic  in  the  reference  wavefield.  We 
decompose  initial  data  into  plane  waves  as  is  suggested  by  Ludwig’s  construction. 
Then  geometrical  optics  is  used  to  propagate  each  plane  wave  through  the  pertur¬ 
bation.  The  plane  waves  are  put  back  together  again,  and  caustics  are  located  by 
the  failure  of  a  stationary  phase  calculation  to  reduce  the  number  of  integration 


variables. 
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3.1  _An  Appropriate  Problem 

To  begin  this  investigation,  we  analyze  the  simplified  problem  of  a  strictly  convex 
wavefront  moving  into  a  two-dimensional  region  where  the  reference  wavespeed  is 
constant,  Co  =  1.  Since  the  wavefront  is  convex,  a  cusped  caustic  will  develop  even 
though  the  region  is  homogeneous.  We  take  8c  to  be  a  localized  high  frequency 
perturbation  in  some  direction  k, 

8c  =  x(x)e,k-x 

where  x(x)  iS  an  envelope  function.  We  will  later  discuss  what  happens  when 
solutions  are  summed  over  k.  See  Figure  3.1  for  an  illustration. 

Next,  we  choose  the  form  of  the  singularity  in  the  initial  data.  The  fundamental 
solution  of  the  wave  equation  in  2-D  with  U  =  0  for  t  <  0  and  c  =  1  is 

which  suggests  using  the  following  initial  data: 

Wo(x,  0)  =  a(x)-^=H(-ip(x)) 

vw) 

where  i/>(x)  =  0  describes  the  initial  wavefront,  and  a(x)  is  smooth,  with  compact 
support.  The  initial  value  ^-(x,  0)  is  determined  by  specifying  that  it  be  a  single 
progressing  wave,  to  leading  order. 


To  be  assured  that  a  simple  caustic  will  occur,  we  assume  that  the  incident 
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wavefront  is  strictly  convex.  To  simplify  calculations,  we  also  assume  that  ^(x)  can 
be  described  by  its  graph,  0(x)  =  x?  —  T(a:i). 


3.2  Decompose  Initial  Data  into  Plane  Waves 


Since  the  geometrical  optics  representation  breaks  down  near  a  caustic,  we  decom¬ 
pose  the  initial  data  into  plane  wave  components,  each  of  which  can  be  propagated 
through  the  region  containing  the  perturbation  Sc  without  forming  a  caustic: 


Then 


£/(x,0)  =  J/o(x,0)  =  f  c19uq(x,0’,0)  . 
J\o\=i 


6U(x,t)  =  J  d06u(x,  t\ 9) 


where 


Sutt  -  V2Su  =  2x(x)eikocuo„ 

(9) 

Su  =  Sut  =  0,  for  t  <  0  . 

To  accomplish  the  decomposition  we  write  the  initial  data  as  the  inverse  of  its 


Fourier  transform: 


<7(x,0)  =  /  d9  I™  duueiu$-xF(U(x,  0)) 
Jm=i  Jo 


where 


F(U(x.  0))  =  j^dy 
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is  the  Fourier  transform  of  £/(x,  0).  Making  the  change  of  variables  77  =  (yt,  - rp(y ))  = 
(2/i>r(y1)  -  y2)  and  evaluating  H(i}2)  gives 

F{U(xt  0))=  r  dm 

J— 00  JQ  yV2 

By  writing  a(n,r(Vl)  -  72)  =  and  invoking 

the  assumption  that  a  and  have  compact  support  in  772,  F(U(x,  0))  can  be 
approximated  by 
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Notice-that  to  each  point  on  the  initial  wavefront  t/»(r/)  =  0  is  associated  a  direction 
0,  which  is  the  unit  normal  to  ip(r])  =  0  at  that  point.  That  is, 

*  =  jWjW)  (11) 

since  the  stationary  point  is  given  by  9  •  (1,  T'^))  =  0  and  rj)(rf)  =  rj 2  —  r(r/i). 


Substituting  the  asymptotic  evaluation  of  F(U(x,  0))  into  the  expression  (10) 
for  U (x,  0)  gives 


The  limits  of  integration  of  the  u>  integral  can  be  extended  to  include  the  negative 
u  axis  by  changing  the  sign  of  both  u>  and  9.  Then 


r  du =  8{x-9-  {tjI r(77*))  •  9), 

J- CO 


and  the  initial  data  is  decomposed  into  plane  wavefronts  with  a  ^-function  singu¬ 
larity, 

t/(x,  0)  ~  ||=1  Mow.  ' »  -  «.  *M»  •  «>• 

Let  t(9)  =  F(vi))  ‘  This  is  the  distance  from  the  line  x  •  9  =  0  to  the 
tangent  to  the  curve  ^(x)  =  0  at  (j/j,  r(j/J))  and  can  be  considered  as  a  time  delay 
when  solving  (9).  (See  Figure  3.2.)  To  keep  the  notation  simple  we  will  drop  t{9) 
until  we  obtain  8u(x,t;0),  at  which  point  we  will  shift  t  — ►  t  —  t(9). 
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3.3  -Geometrical  Optics  Approximation 

We  can  now  solve  each  plane  wave  problem.  Clearly,  the  reference  field  uq  is  of  the 
form  uo(x,  t\ 9)  =  a(9)S(t  -x-9).  Then  a  progressing  plane  wave  representation  can 
be  derived  for  the  perturbational  field  Su,  where 

(w~vjSu= 2*(x)e,1,'xJ?“°  (12) 

This  representation  will  be  valid  up  to  the  time  that  caustics  form  in  the  perturba¬ 
tional  field. 

The  calculation  of  Su  is  very  similar  to  a  calculation  done  by  Symes  and  Santosa 
in  [13].  We  guess  that  Su  has  the  expansion 

Su(x,t;9) = 

bi(x,  t)6'(t  -  tf(x))  +  62(x,  t)S(t  -  <j>{x ))  +  63(x,  t)H(t  -  <f>({x)), 

where  <f>(x)  is  the  incident  phase  <f>(x)  =  x  •  9. 

Substitution  of  this  expression  into  (12)  with  |V^|2  =  |0|2  =  1  gives 

(2 bu  +  2Vbt  ■  V<f>  +  byV2<f>  -  2 x(x)eik-xa)S"{t  -  <j>) 


+  (2 b2t  +  2V62  •  +  62V2^  +  but  ~  V26i)<5,(t  —  <f>) 

+  (26*  +  2V63  ■  +  63VV  +  b2tt  -  V\)6(t  -  <f>) 

+  (63„  -  -  <t>)  =  0 
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A-or(x,o(x))  =  kx 


(17) 
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V^V[^r(x,<£(x))]  =  V^.V[k.x] 
V^r-  V^  +  <£r,|V<£|2  =  k  -V(f>. 


(18) 


Note  that  the  initial,  condition  together  with  the  eikonal  equation  implies  Snell’s 
equal  angle  law  of  reflection.  From  (17), 

f  =  ^  +  v*. 

Yrt  Yrt 

That  is,  the  vector  is  the  sum  of  two  unit  vectors  since  |V<£r|2  =  <£2  and 
|V</>|2  =  1.  The  physical  ray  ^  =  -“T1  is  given  by  the  Hamilton- Jacobi  equations. 
Then 

k  dx 

which  is  precisely  Snell’s  law,  as  can  easily  be  seen  from  Figure  3.3. 


Using  (18)  in  (16)  and  ^(x)  =  x  •  9  gives 


b\  = 


X« 


k  •  V<f>  k  •  6 

The  second  condition  in  (13)  gives  the  following  differential  equation  for  62(x,  t) 
on  the  incident  wavefront: 


262,  +  2V62  •  +  62V2<£  = 

3 bUt  +  4 V6lt VV  +  V26,  on  t  =  <f>(x) 
Assuming  that  62(x,  t)  can  be  represented  by 


(19) 
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and  substituting  this  expression  into  (19)  yields 


2 b%,  +  2 b°2V<j>r  •  V<f>  =  3 b\4>2rt  +  Ab\<j>ry<f>r  •  +  b\\V<t>r\2 


as  the  highest  order  term,  0(k).  We  use  the  eikonal  equation  for  the  reflected  phase, 
=  \V<f>%  which  comes  out  of  the  fourth  condition  in  (13)  and  will  be  derived 
later  in  this  section.  Then 

6”  =  2&i<Mx,<«x))  =  2^<Mx,<A(x)). 

k  •  9 

The  third  condition  in  (13)  gives  an  equation  for  &3(x,i)  on  t  =  <f>(x): 


263,  +  2V63  •  V<j>  +  63  VV  =  -6lttt  -  2V6lt,  •  V<j>  -  blttV2<}>  -  V26lt 

+62,,  +  2V&2,  •  +  b2,^2<j>  +  V2&2 

on  t  =  <f>(x) .  (20) 


Substituting  the  expansion 


gMe.^(x,  0 


(21) 


into  (20)  yields  the  following  as  the  highest  order,  0(k2),  term: 


263-Vrt+263-1V^r-Vd  =  -6{<-26^2V^-V^-611|V^r|2 

+b\£t  +  2b\<j>ri¥ <f>r  ■  +  b°\V<j>r\2  . 
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Again-using  <f>] lt  =  |V<£r|2  we  get 

63 1  =  -b\<}>rt  +  b2<f>r,  on  t  =  <f>(x) 

=  ^r#,(M(x)).  (22) 

k  •  0 

Finally,  the  fourth  condition  in  (13)  is  the  wave  equation  for  h(x,  t)  in  the  light 
cone  t  >  <f>(x ): 

b3„  ~  V2&3  =  0  .  (23) 

Using  the  above  expansion  for  &3(x,  t)  in  (23),  we  get  from  the  0(fc3)  and  0 (k2)  terms 
the  eikonal  and  first  transport  equations: 


4t,  -  |v^r|2  =  0 

2  KK  -  2V63_1  •  W,  +  (K  -  V2^r)6j1  =0,  t>  <j>(x)  . 


(24) 

(25) ’ 


The  solution  to  the  eikonal  equation  (24)  with  the  initial  condition  k<f>T{x,  <f>(x))  = 
k  •  x,  where  <p(x)  =  x  •  9,  is  easily  found  to  be 


d>r(x,<)  =  k,- 


2k  •  9, 


x  + 


2k  -e 


Then  the  transport  equation  (25)  reduces  to 


6r,b£-Vb?-V<j>r=0. 


Since 


db~l 

da 


=  tj'g  +  TOJ' •  £,  «' 


■la' 


is  constant  along  rays  described  by 


dx 

da 

dt 

da 


$ 

2k  -9 


Prt  — 


1 

2k  -0 
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or  ft  =  0  -  2k  •  0k. 

On  the  incident  wavefront,  6jl  is  given  by  (22)  so  that 

^3 1  (*,  0  =  77r~~ra(fl)  x  =  x0  -  (2k  •  0k  -  0)(t  -  0  •  x0) . 

4(k-0)3 

Then,  shifting  the  time  variable  i  to  incorporate  the  time  delay  r(0)  from  the  plane 
wave  decomposition,  the  perturbational  field  can  be  represented  as  follows: 


ct  i  a\  'fc((k — (<-t(5))] 

6u(x,  t‘,8)  ~  e  Jk  o  2k  « 


x(x)°WV(i-r(«)-X-«) 


l  &k  •  0 
X(x)a(0) 


<5(<  -  r(0)  -  x  •  0) 


(k-0)2 

4(k  •  Oy 


} 


where  x  =  Xo  —  (2k  •  0k  —  6){t  —  r(0)  —  x o  •  0)  (or  x0  =  (/  +  (2k  •  0k  —  0)  <S>  0)“ 1  (x  + 
(2k  -0k  -0)0). 


3.4  Sum  the  Plane  Waves 


After  summing  over  the  contributions  from  each  plane  wave,  the  perturbational 
field  at  the  surface  xi  ~  0  is 


6U\.  _0  =  f  dOSu  |  - 

l\6 


|«l=i  4(k  ■  0)3 


x(xo)a(0)e,Wr(x>f-r(0)) 


1X2=0 


since  the  incident  wavefront  moves  away  from  x?  =  0. 


In  order  to  evaluate  the  integral  it  is  convenient  to  change  variables.  There  is  a 
natural  change  of  variables  given  by  the  correlation  (11)  between  a  direction  0  and 
a  point  on  the  initial  wavefront.  We  let  z{zx)  =  (*i,r(*,)),  where  tf(z)  =  0.  Then 
=  jv'0(z)|  =  (I  +  r'/2(-i))~?(-r,(s1),  1)  and 

T(°(2 1))  =  -*(-»)  •  %*i)  =  (i  +  r/2(2l))-1/2(-r(z1)  + 

The  curve  ip(z)  ~  0  is  strictly  convex  so  that  the  map  0  i-»  zx  is  invertible.  This 
gives 


+r(--i)  -  +*,(r(Sl)  +  2(k  -  *,[")*,)) . 
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3.5  -Locate  Caustics 


Caustics  occur  where  the  stationary  phase  expansion  of  (26)  is  invalid,  which  is 
where  the  phase  <f>T{z\)  has  degenerate  stationary  points.  Using 

ST  =  2(*J  -  *,[■')  +  (1  +  r!)"'/!ri  _  ** +  *,(1  _ 

=  0 


gives 


(P<f>r 


T" 


r 


[(1  +  r«)*/a 


t-1 


dz\  0(k2  -  kiP) 

Thus  the  caustic  locus  of  the  perturbational  field  on  the  line  x2  =  0  is 

1 


t  = 


r*(*i) 


(i  +  r'2(i,))3/2. 


(27) 


Note  that  it  is  independent  of  the  perturbation  direction  k  and  the  receiver 
location  x\.  Also  notice  that  the  right  hand  side  is  the  radius  of  curvature  of 
ij>(z)  =  0.  This  leads  one  to  think  that  (27)  should  also  represent  the  caustic  locus 
of  the  incident  wavefront.  This  is  verified  by  studying  the  map  (zi,<)  *-*  (zi,^)  by 
which  points  on  the  initial  wavefront  are  mapped  along  rays  to  the  wavefront  at  a 
later  time  t.  That  is, 


(*i.  *,)  =  (--i,r(2,))  +  f(i  +  ra(*i))-,/,(-r'(*.M). 


A  caustic  is  found  by  looking  for  where  this  map  is  singular.  The  Jacobian  is 

det  =  (Hz,)  +  i)-‘((r'2(z,)  + 1)3/2 -«r"(z,)], 
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which-is  zero  precisely  when 

,=4o(l+rW’- 

Thus,  if  there  is  a  caustic  in  the  perturbational  field  at  the  surface  £2  =  0,  then 
the  distance  along  a  ray  from  the  incident  wavefront  to  the  incident  caustic  must 
equal  the  distance  along  a  ray  from  the  incident  wavefront  to  the  perturbation  plus 
the  distance  along  the  corresponding  reflected  ray  to  the  perturbational  caustic. 
This  can  only  happen  when  the  perturbation  6c  =  x(x)e'k‘x  is  located  above  the 
incident  caustic.  (See  Figure  3.4.)  In  that  case,  it  is  possible  to  represent  the 
scattered  field  at  the  surface  using  the  methods  of  Ludwig. 

If  the  perturbation  is  located  near  the  incident  caustic,  the  only  caustic  present 
in  the  perturbational  field  is  that  piece  of  the  incident  caustic.  Instead  of  amplifying 
the  perturbational  field,  the  caustic  causes  the  reflected  rays  to  be  spread  out.  (See 
Figure  3.5.)  If  the  perturbation  is  below  the  incident  caustic,  there  will  be  no  caustic 
at  all  in  the  perturbational  field.  In  fact,  the  incident  wavefront  is  concave  there, 
so  that  the  reflected  wavefront  is  also  concave.  (See  Figure  3.6.)  For  the  cases 
above,  the  scattered  field  at  the  surface  can  be  represented  by  a  geometrical  optics 
expansion.  In  any  case,  8U  can  be  represented  using  known  methods,  even  if  the 
perturbation  sits  on  a  caustic. 
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*3  =  0 


initial  wavefront 


caustic 


Figure  3.5.  The  perturbation  is  located  on  the  caustic  in  the  reference  field. 
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3.6  -Integral  Operator  Representation  for  the 
Scattered  Signal 

If  there  are  no  caustics  in  the  perturbational  field,  then  the  method  of  stationary 
phase  can  be  used  to  asymptotically  represent  the  integral  in  SU  (26).  Without  loss 
of  generality,  we  can  assume  that  there  is  only  one  stationary  point.  We  then  get 
the  following  form  for  the  scattered  signal: 

SU |«=o  ~  , 


where 


and  <f>r(k,xi,t)  =  k4>T{z\,x\,t, k). 


X?  =0 


To  generalize  this  formulation  to  general  perturbations,  Sc,  we  can  take 


Sc  =  x(x)£c(k)e‘ 


ik-x 


to  be  one  localized  Fourier  component  of  a  general  velocity  perturbation  5c(x). 
Then  the  previous  analysis  will  still  be  valid. 

Since  Sc  occurs  only  in  the  right  hand  side  of  the  first  transport  equation, 

2bu  +  2V&i  •  v<j>  +  bi  v2(j>  =  2x(x)dc(k)elk'xa 


we  can  make  a  simple  modification  to  account  for  more  general  perturbations.  We 
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need  only  change  the  ansatz  for  &i  to  account  for  a  factor  of  £c(k): 

n=l  \lK) 

The  leading  order  term  in  61  is 


6}  =  -£Lfc(k) . 
1  4k -0  K  f 


Similar  expressions  can  be  obtained  for  bi  and  63.  The  expression  for  the  scattered 
field  will  change  only  by  a  factor  of  5c(k): 


8U\S  =0  ~  /l(k,  xu  <)e'^r(k’ri,t)5c(k) 


Summing  over  k  gives  the  trace  of  the  perturbational  field  8W  due  to  a  localized 
velocity  perturbation  8c(z), 

8W\X1=0  ~  J dkA(k,xltt)ei^k'Xl>t^(k) . 

We  can  write  £c(k)  as  a  Fourier  integral, 

W|XJ=0  ~  f  dk  A(K  *1,  0e^r(k,x,’°  J  dz  8c{ z)e-'k-z  . 

Then  the  scattered  signal  can  be  represented  by  the  following  Fourier  integral  op¬ 
erator  acting  on  the  velocity  perturbation  8c( z): 


[pfc](i„«)= 

j  J  <fk<fzA(k,  i,,  ()e'|,f'<kt'"‘1“k'z|Sc(z) . 


CHAPTER  4 


An  Example  of  Anomalous  Scattering  Strength  from  a  Caustic 

In  the  previous  chapter  it  was  determined  that  if  a  high  frequency,  unidirectional 
perturbation  were  located  near  a  caustic  in  the  incident  field,  the  energy  would  be 
spread  out  instead  of  focused.  This  result  does  not  extend  to  general  perturbations, 
Sc.  It  was  expected  that  the  product  of  Sc  and  uo„  in  the  right  hand  side  of  the 
perturbational  problem 

Sutt  —  V26u  =  2  Sc  i/o„ 

would  cause  the  scattered  signal  to  be  illuminated  if  Sc  were  located  near  a  caustic 
in  the  reference  field  «o-  This  is,  in  fact,  the  case.  Caustics  in  a  probing  wavefield 
do  cause  velocity  perturbations  to  scatter  signals  more  strongly  than  a  smooth  ray 
pattern  would. 

In  order  to  quantify  this  statement,  we  find  a  sequence  of  general  perturbations 
Sc  that  generate  0(  1)  scattered  signals  in  a  smooth  reference  field  and  increasingly 
larger  signals  in  the  presence  of  a  caustic.  That  is,  — ►  oo  as  Sc  moves  closer 

and  closer  to  the  caustic.  This  example  illustrates  the  fact  that  even  if  there  is  no 
caustic  in  the  perturbational  field,  the  scattered  signal  can  be  anomalously  large. 
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The  construction  of  this  example  is  best  described  by  the  pictures  in  Figure  4.1. 
A  perturbation  Sc  is  constructed  so  that  the  incident  wavefront  will  be  reflected  as 
a  plane  wave.  The  perturbation  is  localized  to  lie  within  a  ray  tube  in  the  incident 
field.  As  Sc  moves  closer  to  the  caustic,  the  reflected  ray  tube  gets  smaller.  That 
is,  the  same  amount  of  energy  is  reflected  into  a  decreasing  area. 

We  first  construct  the  scattered  field  due  to  a  velocity  perturbation  of  the  form 

Sc  =  X{x,z)eik'l{x’z) . 

The  incident  field  at  t  =  0  is  taken  to  be 

uo(x,2,0)  =  a{x,z)H(-rp{x,z)) , 
which  gives  the  reference  field 

uo(x,  t )  ~  a(x,  z,  t)H(t  —  <f>(x ,  z)) . 

The  phase  <p  satisfies  the  eikonal  equation  |  V<£p  =  1,  and  the  amplitude  a(x,z,  t) 
is  given  by  the  transport  equation 

2Va-  W  +  aVV  =  0. 


The  Heaviside  singularity  H{-ip)  is  chosen  rather  than  the  point  source  singularity 
-^=H(—Tp)  because  the  construction  is  simpler.  In  fact,  it  parallels  the  construction 
in  Chapter  3. 


perturbational  wavefront 


Figure  4.1.  As  the  perturbation  is  placed  closer  and  closer  to  the  caustic, 
the  reflected  ray  tube  gets  smaller  and  smaller. 


perturbation*!  wavefront 


caustic 


initial  wavefront 


wavefront 


Figure  4.1.  As  the  perturbation  is  placed  closer  and  closer  to  the  caustic, 
the  reflected  ray  tube  gets  smaller  and  smaller. 


perturbational  wavefront 


Figure  4.1. 


As  the  perturbation  is  placed  closer  and  closer  to  the  caustic, 
the  reflected  ray  tube  gets  smaller  and  smaller. 
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The  scattered  field  has  the  form 


where 


6u(x ,  2,  t )  =  b\(x,  z,  t)8(t  —  <f>(x,  2))  +  62(x,  2,  t)H(t  —  <f>(x,  2)) 


The  phase  must  satisfy  the  initial  condition 


Or(l,2,^(x,2))  =  7(X,2)  (28) 

and  the  eikonal  equation 

|V^r|2  =  <  • 

We  want  to  construct  the  perturbations  5c  =  x(*>  2)e,fc7^’^  in  such  a  way  that 
the  reflected  phase  <f>T  is  a  plane  wave  that  moves  in  the  —2  direction.  That  is,  we 
want 

$r(®j  “>  0  —  t  Z  * 

This  automatically  satisfies  the  eikonal  equation.  In  order  to  satisfy  the  initial 
condition  (28),  we  must  take 

7(1,2)  =  <j>(x,z)  +  z  . 

The  coefficient  b\  satisfies  the  transport  equation 

2b\ +  26}  V<pT  ■  V<f>  =  2 xaeik{y-*r) .  (29) 
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Taking  the  gradient  of  the  initial  condition  (28)  with  -y(x,  z )  =  <f>(x,  z)  +  2  gives 
V<f>r(x,  Z,  <f>(x ,  z))  +  <f>Tt{x,  2,  <j>(x,  z))V<f>{x,  z)  =  V<£(x,  z)  +  (0,  1) 

and 

V<£r  •  V<£  +  <f>u \V<f>\2  =  \V<f>\2  +  <i>t . 

But  |V^|2  =  1  so  that 

V4>r  ■  V<f>  +  <f>r,  =  1  +  <t>* .  (30) 

Using  (30)  in  (29),  we  obtain 

,1  *a  Af  \  * 

Oi  =  — 7-  on  <f>(x,z)=t. 

1  +<Pz 

Similarly, 

b°2  =  2&|<£r,(x,2,<£(x,2))  =  r^jr-  on  <j>(x,z)  =  t. 

i  +  <Pz 

The  coefficient  62  satisfies  the  wave  equation  in  the  light  cone 

62,,  —  V262  =  0  in  t><f>(x,z).  (31) 

Substituting  the  above  expansion  for  62  in  (31)  gives  the  eikonal  equation  for  the 
reflected  phase  and  the  following  transport  equation  for  6°: 

^,6°,  -  •  Vd>r  =  0  in  t><}>(x,z). 

Then  6°  *s  constant  along  the  rays  described  by 

I  = 

<tt 
da 


-K  =  -1 
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2=0 


or  jfj=  (0,-1),  where  x  =  ( x,z ).  Then  the  reflected  field  at  the  surface  can  be 
written  as 

l+P*(*0,2o) 

where  the  rays  are  given  by 

(x,z)  =  (xQ,</>(x0,z0)  -  t)  , 

and  (xo,  z0 )  is  on  the  incident  wavefront  <f>(x o,  zq )  =  c.  Notice  that  the  amplitude 
of  the  reflected  field  is  constant  along  reflected  rays.  This  construction  is  valid  only 
above  the  incident  caustic,  but  this  is  sufficient  for  present  purposes. 

To  form  the  sequence,  let  each  element  of  {5cn}  be  a  localized  high  frequency 

* 

perturbation  of  the  form 

f>Cn  =  Xn('*»  z)e,fc7(x’*')  , 


where 


7  (*,*)  =  4>(x,z)  +  z 


In  order  to  describe  the  characteristic  function  Xn(^i~)i  let  T  be  a  slim  ray  tube 
enclosed  by  two  incident  rays.  The  cross  sectional  area  of  the  ray  tube  goes  to  zero 
as  the  rays  approach  the  caustic  (see  Figure  4.2).  Let  r  be  the  length  of  the  bisector 
of  the  two  bounding  rays  between  the  initial  wavefront  and  the  crossing  of  the  rays, 
as  in  Figure  4.3.  Then  let  \^(x.z)  be  equal  to  one  on  a  disk  T>n  that  is  centered 
at  a  point  on  the  bisector  a  distance  r  —  from  the  incident  wavefront.  This  disk 


Figure  4.2.  T  is  a  slim  ray  tube  enclosed  by  two  incident  rays. 


Figure  4.3.  The  length  of  the  bisector  of  the  two  bounding  rays  between  the 
initial  wavefront  and  the  crossing  of  the  rays  is  r. 


Figure  4.4.  The  envelope  function  x„(x)  is  one  on  the  disk  T>n. 
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should- be  contained  strictly  within  the  ray  tube,  and  Xn(^)  should  be  zero  outside 
T  (see  Figure  4.4).  As  n  gets  larger,  Vn  moves  closer  to  the  caustic. 

The  relative  size  of  the  scattered  signal  to  the  perturbation  is  given  in  terms  of 
the  L 2  norms: 


||«C„||2  =  J  J  dxodzo\Xn(x<h  Zo)\2  =  \\Xn\\2 
=  JJdxdtX„(x  0,2o)ri£gL5 


z=0 


(32) 


The  time  that  it  takes  to  get  from  the  initial  wavefront  to  the  reflector  then  to  the 
surface  is  <^(xo,2o)  +  ~o  so  that  points  (x0>2o)  are  mapped  to  (x,  t)  by 

(x,  t)  =  (xo,<f>(x0,  z0)  +  z0)  at  z  =  0. 

This  is  illustrated  in  Figure  4.5.  Making  this  change  of  variables  in  (32)  changes 
the  amplitude  only  by  the  Jacobian  factor  J  =  (1  +  <bx).  Then 

a(x0,  -s'o) 


II  || 2  =  J  J  dxodzo 


Xn(Zo,  Zo) 


\J  1  +  <&r(ZO>~o) 


XnO. 


\/T+lfz 


Since  <j>  satisfies  the  eikonal  equation  |V<£|2  =  1  and  <f>x  >  0,  the  factor  (l  +  <£*)_2 
is  0(1).  But  as  the  incident  rays  get  closer  to  the  caustic,  the  amplitude  a(x0,  zo) 
approaches  infinity.  This  was  shown  in  Chapter  2  by  integrating  the  transport 
equation  over  the  ray  tube,  which  collapses.  Thus 

w 

\\Scn  |P 


— >  oo  as  n  — >  co  . 


Figure  4.5.  Points  (x0,  zo)  are  mapped  along  reflected  rays  to  (x\  t)  at 
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If  there  is  no  caustic  in  the  incident  wavefield  then  is  always  bounded  since 
the  amplitude  remains  bounded. 

An  incident  point  source  field  could  have  been  used  instead  of  the 

Heaviside  field,  H(—ip).  Under  these  circumstances,  the  conclusion  would  be  that 

iMhi 


Pell 


remains  bounded 


(33) 


as  long  as  the  support  of  8c  stays  away  from  the  incident  caustic,  but  approaches 
infinity  if  the  support  of  Sc  approaches  the  caustic. 

This  is  verified  by  a  result  of  Rakesh  [16],  which  states  that  if  there  are  no 
caustics  present,  then 


iK-U-u*  -  constantll,5cll' 


(34) 


The  s-norm,  ||  •  ||s,  is  the  Sobolev  norm  defined  by 


ii/ii.  =  (/Ki><fk(i+t2)'i/(k)i: 


In  two  dimensions,  (34)  can  be  written  as 

S*«U-Li  -  constantl^cll,  , 

which  says  that  the  s  -  j  th  derivative  of  6u  is  the  same  size  as  the  sth  derivative 
of  6c.  Since  the  factor  of  lc~'i  in  (33)  acts  like  a  (—5  )th  derivative,  the  results  (33) 
and  (34)  are  consistent. 


CHAPTER  5 


The  Inverse  Operator 


5.1  The  Construction  of  the  Inverse 


In  Chapter  3  we  found  the  following  Fourier  integral  operator  representation  for  the 
scattered  signal  due  to  a  general  perturbation  6c( z): 


(35) 


*n2=0  ~  (p*ci(*i»o 

=  J  J dkdzA(ktxltt)e^^^^Sc(z) . 

This  representation  is  valid  away  from  any  caustics  in  the  perturbational  field.  It 
is,  in  fact,  valid  even  if  there  are  caustics  in  the  incident  field.  The  goal  of  this 
chapter  is  to  find  an  inverse  operator  for  this  representation.  That  is,  we  will  find 
an  inverse  of  the  forward  map 

£>S(c„|«c  =  SW\„^  ■ 


Let  Pi  be  the  operator  defined  by 


PM*  i,0  =fj  dk-x  dz.4i  (kj ,  1 1 ,  t  )ci^,(kl  ,Xl  ’zl5c(z) 


The  operators  P  and  Pi  are  the  same;  the  subscript  1  is  used  only  for  clarification. 
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We  intend  to  find  an  operator  P2  that  is  close  to  an  inverse  operator.  To  define 
what  close  means,  we  must  introduce  a  parameter  r  that  describes  the  roughness 
of  Sc.  Recall  that  we  assumed  in  the  geometrical  optics  calculations  that  Sc  acted 
as  a  reflector.  This  is  a  statement  that  Sc  is  rough.  One  way  of  viewing  this  is  if 
the  velocity  perturbation  can  be  written  in  the  form 

Sc( z)  =■  a(z)e‘T'1'^  . 

As  r  — +  00,  Sc  oscillates  more  and  more  rapidly  and  thus  becomes  very  rough. 
Assuming  that  Sc{z)  can  be  represented  in  this  way  and  a(z)  has  compact  support 
in  z,  we  will  look  for  an  operator  P2  such  that 

P2P\Sc{z)  ~  ISc{z)  as  t  — ►  00  . 

Define  P2  as 

P2/(y)  =  JJ  dk2dxA2(k2,.x)e-^kj’x)-kj-yJ/(x) .  (36) 

P2  is  similar  to  the  adjoint  operator  to  P\.  Then,  letting  x  represent  (x\ ,t), 

P3PlSc( y)  =  1 1 J  J  dkidkidxdze'M*' .*) -ki  .*)+k* -yl 

x  Ai(ki,x)A2(k2,x)£c(z) . 

Since  Sc(y)  can  be  written  as  the  inverse  of  its  Fourier  transform, 

Sc(y)  =  J  J  dkdzelk'^y~z^Sc(z)  , 


(37) 


59 


the  idea  is  to  make  PiP\8c( y)  look  like  (37).  In  order  to  accomplish  this,  two  of 
the  integrals  in  P2Pi#c(y)  must  be  eliminated.  But  a  large  parameter  in  the  phase 
function  is  necessary  to  perform  a  stationary  phase  calculation.  Therefore,  we  make 
the  change  of  variables 

k,-  -*  rk,-,  i  =  1,2  . 

The  amplitude  in  A\  was  shown  in  Chapter  3  to  be  homogeneous  of  order  £  .  Then 
assuming  that  Ai  is  homogeneous  of  order  £,  we  have 

P2Pi<Sc(y)  =  J  J  J 

x.Ai(ki,x)A2(k2,x)6c(z) , 

where  $(ki,k2,x,y,z)  =  <£(ki,x)  -ki  •  z  -  <£(k2,x)  +k2  -y. 

We  intend  to  do  a  stationary  phase  calculation  on  the  (k2,x)  variables.  There¬ 
fore,  the  error  term  will  be  expressed  as  an  integral  over  the  (kj ,  z)  variables  of  a 
function  that  is  0(r-1).  In  order  to  ensure  that  the  entire  integral  is  0{ r-1)  and 
can  be  considered  asymptotically  of  lower  order  than  the  leading  term,  we  do  the 
following  calculation.  The  integral  P2PiSc( y)  can  be  rewritten  as 

PiPMy )  =  J  J  dkidzr4+*+<a(k,,z)e,rk'  za(z)eiT7(z) 

where 


a(ki,z)  =  J  J  t/k2dx/li(kt ,  x)/l2(k2, x)e*r^(k,,i,_^k2'^+k2’y^ 
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Then,-letting  6(ki,z)  =  a(z)a(ki,z)  and  q  =  4  +  |  +  £, 

P2PM  y)  =  J  j  dk1dzr?6(k1,z)e‘'Tkl-ze,'T7(z) 

=  f  [  dk1dzr76(k1,z)e,vk>-ze,T1,(z) 

*/  Jtj  <  1 

+  J  Jk  >idk1dzr’6(k1,z)e,Tkl'ze,'T7(z) . 

The  amplitude  6(ki,z)  has  compact  support  in  z  since  a(z)  does.  Therefore,  the 
integral  over  the  region  ki  <  1  is  bounded  by  something  that  is  0(r7), 

/  /  dk1dzr’?6(k1,  z)e'Tkl-ze,T7^z^  <  constant  O(r’)  . 

J  J k\  <1 

For  the  infinite  integral,  we  write 


irk|  Z  _  V7  -irki  Z 

e  - 


then  do  an  integration  by  parts  in  z: 

/  [  dkldzr96(k1,z)e,r7(z)T^7-V.-e,'Tkl  Z 

=  -  f  f  dMzr9^  •  [v.7(z)6(k1,z)  +  J-V^k^z)]  e‘^<z)e‘>kl  Z  . 

J  Jki>  1  Ki  L  IT  J 

Thus,  one  integration  by  parts  pulls  in  a  factor  of'fcf 1  without  changing  the  order  in 

r.  Two  more  applications  of  integration  by  parts  provides  a  factor  of  k f3  ensuring 

the  absolute  convergence  of  the  integral,  again  without  changing  the  order  in  r. 

It  is  possible  now  to  estimate  the  integral  PiPi8c( y)  as 

P2PM y  <  Tq  J  J  dkjdz  |&(ki,z)e’rk,’ze,nf^| 

=  0(r<). 
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Now  stationary  phase  can  be  performed  on  the  (k2,x)  variables.  The  stationary 
points  are  given  by 

=  y- Vfcj^(k2,x)  =  0  (38) 

=  V^(kltx)-V^(k2,x)  =  0.  (39) 

The  second  equation  (39)  is  satisfied  by  k2  =  ki.  The  first  equation  (38)  specifies 
some  stationary  points  (kj,xJ).  As  long  as  k2  =  ki  is  the  only  solution  of  (39)  then 
the  following  calculation  of  the  approximate  inverse  can  be  done.  It  doesn’t  matter 
how  many  points  xj  there  are.  For  simplicity  of  notation  we  will  assume  that  there 
is  only  one;  call  it  x*.  In  any  case,  the  implicit  function  theorem  ensures  that  the 
stationary  points  will  be  locally  unique  as  long  as 

det|VfclV^(kt,x*)|^0. 


This  mixed  Hessian 


H  = 


firiki  xikt 


is  singular  precisely  at  a  caustic  and  corresponds  to  the  singularity  of  Beylkin’s  ma¬ 
trix  (8).  If  H  is  non-singular,  the  order  in  k  of  the  amplitude  Ai(ki,  x',  t)  determines 
the  size  of  the  operator  Pi  (15). 


Notice  that  condition  CM)  concerns  the  ray  direction  vector 


/ dii  dx2  dt  \ 

l  d7’  d7’d^) 


do  <)o  dd 
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Thus,_the  assumption  that  the  solution  of  (39)  is  globally  unique  is  a  statement  that 
the  direction  of  the  rays  at  the  surface  change  when  the  direction  of  the  perturbation 
changes.  For  an  illustration  when  the  background  velocity  is  constant,  see  Figure 
5.1. 

At  the  point  (k2,x)  =  (ki,x‘),  the  second  derivatives  of  $  are  given  by 

VjtjV*j$|(k,,x*)  =  -VjfcjVfcj^ki.x*) 

ViVi$(kli*.)  =  VfVf^(ki,x*)  -  VfV£^(ki,x,‘)  =  0 
^^2 =  —  X  )  . 

Therefore  det[V(Jtli£)V(^)£)$(k1,x*)]  =  -det[V*,  V*^(ki,x”)],  and 

P*PMy)  ~  J  j dk1dzr2+?+r(2T)2Ai(ki,x*)i42(ki)x,')(Jc(z) 

e~'i  .  ,  v. 

X  e»r(y-z)-ki 

\/detvj^V£^(ki^k*j 

+  J  J  dk1dze,r^y-z^k'0(r2+2 +<_1)  , 

where  the  A  j  are  eigenvalues  of  V(*2i£)V(jt2i£)$(ki,x*). 

Lettingl  A2(k2,x)  =  ^  [A1(k2,x)]“l  e{f  E^^/det  V*2V£^(k2,x)  and 
£  =  —  i  ,  we  have 

PtPMy)  =  // dk1dzr26c(z)e,'T(y-z)-k‘ 

+  J  J  dkidzSc(z)eiTiy-z)‘klO(r)  . 


Figure  5.1.  The  direction  of  the  reflected  ray  at  the  surface  changes  if  the 
direction  of  the  perturbation  changes. 
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Changing  variables  again 


rk,-  — ►  k; 


gives 

P2PiSc( y)  =  J  J  dkidz  <5c(z)e‘(y~z)'k'  +JJ  dk.dz  Sc{z)e^-z)^0{T~l) . 
As  we  saw  above,  the  error  term  is  therefore 


PiPiSc(y)  ~  /£c(y) . 


5.2  An  Extension  for  3-D  Varying  Background  Velocity 

Note  that  the  representation  of  the  trace  of  the  perturbational  field  SU  found  in 
Section  5.1  (35)  is  general  enough  to  include  the  case  of  a  general  initial  wavefront 
moving  into  a  medium  with  a  smoothly  varying  background  velocity  c0  =  co(x).  In 
this  case  the  reference  field  has  the  form 


u0  =  a(x,  t)8(t  -  ^(x)) . 


Then  the  reflected  phase  <pr(k,x,t)  satisfies  the  eikonal  equation 


5 frA  -  'v^i2  = 0  ■ 


6r(x,4>(x))  =  k  •  x  , 


and  the  amplitude  A  solves  the  transport  equation 


2 4>r,At  -  2 V A  ■  V6r  +  (&„  -  VVrM  =  0.  t>  t(x) 
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To  extend  the  calculation  to  three  dimensions,  note  that  the  vector  x  can  be 
n-dimensional 

x  =  (it, 22, •••,£„),  n  >  1  . 

Let  x  =  (*i,  X2,  •  •  • ,  £n-i).  Then  the  forward  map  can  be  expressed  as 

6W\Xn=o  ~  PMx\t)  =  J  J dkidzA1(k1,x',t)e,^kl'x''t)-k‘'z)^c(z) . 


The  amplitude  Ax  is  homogeneous  of  order  m,  where  m  depends  on  the  dimension 

n  and  the  initial  singularity.  Furthermore,  the  inverse  operator  P 2  can  be  defined 
as 

P2  /(y)  =  J  J  dk2dxA7{k2,x)e-i^k^'^f(x) , 

where  A2  is  defined  as  in  Section  5.1. 


CHAPTER  6 


The  Relationship  to  Beylkin’s  Inverse 

Although  the  inverse  operator  (36)  derived  in  Chapter  5  appears  to  be  different 
from  Beylkin’s  inverse  (8)  described  in  Chapter  2,  they  should  be  related  if  the 
initial  wavefronts  are  the  same.  The  differential  equation  that  Beylkin  treats  is  the 
reduced  wave  equation,  which  is  related  to  the  time  dependent  wave  equation  by 
the  Fourier  transform.  In  this  chapter  we  show  that  the  two  forward  operators  are 
equivalent  up  to  amplitude.  It  then  follows  that  the  inverses  are  also  equivalent. 

The  integral  operator  representation  for  the  scattered  signal  generalized  to  n 
dimensions  is 

[Pi<5c](x',  t)  =  J  J dki  dzAi(k\,x\t)e'k,^r^l,x'’^~^r^6c(z) .  (40} 

The  scattered  field  solves  the  perturbational  wave  equation 

1  a  26c 

—butt  -  V  bu  =  — 3-^0,* 

Co  Co 

bu  =  but  =  0  for  t  <  0, 

where  the  incident  wave  field  u0  satisfies 

1  ml 

~2U0„  -  V  Uq  =  0 
c0 
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uo(x,0)  =  g(ip(x)) , 
and  g(if>(x))  is  singular  at  0(x}  =  0. 

The  Fourier  transform  of  the  integral  operator  representation  for  the  scattered 
signal  found  in  Beylkin  [13]  is 

Fv,c{£,  t)  =  j  J  dxdw  j£^/(x)a(x,  ^)e^(^)+^x.«)-0  .  (41) 

The  scattered  signal  v3C  solves  the  reduced  wave  equation 

u2nlv3C  +  V2u5C  =  — a;2/(x)um 
v3C  =  0  ,  x  large  , 

where  v'n  satisfies 

w2nlvm  +  V2^"  =  S(x  —  tj) 
v,n  =  0  ,  x  large  . 

The  scattered  field  v3C  is  solved  by  the  use  of  the  Green’s  function 

aj2n20vout  +  V  Vut  =  S(x  -  0  . 

If  the  initial  wavefront  u0(x,  0)  were  due  to  a  point  source,  then  the  integral 
operator  representations  of  the  scattered  field,  (40)  and  (41),  should  be  very  closely 
related.  In  order  to  exhibit  this  relationship,  we  first  analyze  the  representation 
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(40)  shown  in  Chapter  5.  We  can  split  the  integral  over  ki  into  an  integral  over  the 

A 

unit  vector  ki  and  an  integral  over  the  length  k\.  Then 

Pl8c(x',t)  =  I J dkrfzAxikux' ,t)8c(z)  J dk, 

if  i4i(ki,x;,  t)  is  homogeneous  of  order  m. 

Since  the  Fourier  transform  of  the  jth  derivative  of  6(x)  is  —(—iky, 


PiSc(x',t)  = 

-J  J  dk1dzi1(kl.x\0<5c(z)(-01'm-n6(m+n-1)(^r(kx,x,,<)-k1-z)  , 


where  <5(m+n_1)  is  the  (m  n  —  l)t/l  derivative  of  8(x). 

This  representation  of  the  scattered  signal  is  an  integral  of  the  form 

1  =  J  dx  a(x)8{m+n-1\<f>(x)) . 

Changing  variables  to 

y  =  (i/l>  •  *  •  >  Vn)  ~  (%1 1  •  •  •  >  %n- 1)  <^(x)) 


gives  the  Jacobian  factor 


J  = 


H  ,  v  ,  u 


-1 


The  solution  of 


>Jn  -  <?(yi,....,yn-uXn)  =  0 
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gives  a?n  =  xn(yn  yn-i)  since  —  is  never  zero.  Therefore 


I  —  [  »„  q(yi»-;-.yn-i,gn(y))  r(m+„_i)/  \ 

j  y §£;{yu---,yn-i,Xn(y))  {Vn) 


which  is  reduced  to 


I  = 


(U 

m+n-1 

i 

*  i  2/n-l?  ^n(y)) 

\dy»j 

a*  (,. 

•  •  >  J/n-1)  ^n(y)) 

=  I  dx  w(x)  . 

J 

Thus  the  scattered  signal  can  be  written  as 

P\6c(x',t)=  if.  .  dkidzu>(ki,x',t,z) . 

J  J<t>r( ki  tX',f)=k| ’Z 


The  reflected  phase  ^r(k|,  x,  t )  is  constant  along  rays,  and  if  there  are  no  caustics 
present,  then  the  rays  do  not  cross.  Thus,  given  (ki,x',f),  there  is  a  unique  point 
z  for  which 

<£r(ki,x\i)  =  ki  •  z 

holds.  The  z  that  satisfies  this  is  specified  by  the  initial  condition  for  the  reflected 
phase 

M ki,z,«^(z))  =  ki  z 


and  is  the  reflection  point  of  the  ray  which  passes  through  (x',  t)  (see  Figure  6.1). 
Since  this  is  an  equation  in  ki  and  z  only,  it  specifies  ki  as  a  function  of  z.  Thus 
(42)  is  an  integral  over  the  normal  vectors  kj,  which  parameterize  the  constant 


z 


Figure  6.1. 


The  point  z  is  the  reflection  point  of  the  ray  which  passes 
through  (x',t). 
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traveL-time  ellipse  U  +  t2  =  constant  (see  Figure  6.2).  If  the  reference  velocity  is 
constant,  then  the  curve  t\  +  £2  =  constant  is  an  ellipse,  otherwise  it  is  ellipse-like. 

Beylkin’s  representation  for  the  scattered  field  can  similarly  be  written  as  an 
integral  over  a  constant  phase  surface, 


-  L 


The  phase  <f>  satisfies  the  eikonal  equation 

I  W|2  =  nl 

which  is  solved  by  the  Hamilton- Jacobi  equations 


dxw(x,^rj) . 


x  =  p 


p  =  -2rc0Vn0 


<t>  =  nl 


Then 

<f>  —  j  TIq  , 

Jrays 

where  the  rays  are  given  by  x  =  p.  Thus  <^(x,  77)  is  the  geodesic  distance  along 
the  incident  ray  from  the  source  point  77  on  the  surface  to  the  interior  point  x, 
and  <£(x,  £)  is  the  geodesic  distance  along  the  reflected  ray  from  the  point  x  to  the 
receiver  position  if  (see  Figure  6.3).  Therefore,  v3C  is  an  integral  over  the  same 
constant  travel  time  ellipse  as  6u.  That  is,  up  to  amplitudes,  the  solutions  to  the 
forward  problems  6u  and  v,c  are  equivalent.  It  is  expected  that  the  amplitudes  will 
also  prove  to  be  equivalent. 
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(*'><) 


Figure  6.2.  The  vectors  ki  parameterize  the  constant  trr*'el  time  ellipse 
tj  + =  constant. 


Since  the  inverse  operators  for  both  Su  and  v3C  are  closely  related  to  the  adjoint 
operators,  the  equivalence  relation  found  for  the  forward  maps  extends  to  the  in¬ 
verses.  The  inverse  operators  for  both  Su  and  v,c  are  also  integrals  over  constant 
phase  surfaces.  In  this  case,  they  are  constant  depth  hyperbolas  with  integration 
over  arrival  times  (see  Figure  6.4). 
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